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Abstract—Results of measurements of the bidirectional reflec- 
tion function of snow for the solar zenith angle close to 54° are 
compared with a recently developed snow optical model based on 
the representation of snow grains as fractal particles. The model 
has a high accuracy out of the principal plane for the observation 
zenith angles smaller than 60°. However, the accuracy is reduced in 
the principal plane. Specular light reflection by partially oriented 
snow plates on the snow surface not accounted for by the model can 
play a role for measurements in the principal plane. The model dis- 
cussed can be used for the grain size retrieval using both ground 
and spaceborne measurements of the snow reflectance. This is sup- 
ported by a high accuracy of the model in a broad spectral range 
545-2120 nm as demonstrated in this work. 


Index Terms—Light reflectance, light scattering, radiative 
transfer, remote sensing, snow. 


I. INTRODUCTION 


OLAR areas are covered by ice and snow. Also, snow 

covers great portions of Eurasia and Nothern America 
in winter months, changing the planetary albedo. It has been 
argued that the snow albedo is modified by an anthropogenic 
interference [4]. 

In particular, the average concentration of soot in the Arctic 
has been increased almost eight times during the time period 
1983-1998 [4]. Such findings are of a great importance and 
should be confirmed by independent measurements especially 
those involving data from orbiting satellites. 

Aoki et al. [3] showed that the snow grain size d and the con- 
centration c of impurities in snow are very important factors for 
the snow albedo r. Both c and d have been retrieved from satel- 
lite data using the algorithm in which snow grains are assumed 
to be spherical [5]. 

To decode the information contained in the spectral signatures 
of snow as detected from a satellite, a comprehensive theory of 
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snow optical properties is needed. This is an extremely complex 
task due to the complex internal structure of snow [13], [16]. 
The horizontal and vertical inhomogeneity of snow (e.g., snow 
layers deposited at different times) add even more complexity 
to the problem. 

The main task of this paper is to validate an approximate 
asymptotic theory (AAT) of snow optical properties developed 
elsewhere [10] using in situ measurements of snow reflectance 
in visible and near-infrared. In particular, measurements at the 
wavelengths 545, 1050, 1240, and 2210 nm have been com- 
pared with the AAT. They coincide with atmospheric windows, 
where gaseous absorption can be neglected or easily taken into 
account. Correspondent measurements have been performed in 
Hokkaido, Japan on February 9, 2001. The AAT is valid only 
for optically thick weakly absorbing snow layers. Therefore, the 
case of extremely dirty snow cannot be considered in the frame- 
work of the AAT. 

The results of this work can be used for monitoring of snow- 
pack properties (e.g., the ice grain size and the concentration of 
pollutants) using satellite measurements. 


Il. THEORY 


The snow reflection function can be modeled using the fol- 
lowing approximate analytical solution of this equation valid in 
the limit of small light absorption by snow grains [10] 


RQ) = Ro(3) exp (~a f (0) . a) 


Here R is the bidirectional reflectance (or the reflection func- 

tion) of a semi-infinite snow layer. The vector-parameter @ has 

coordinates 9o, 0, y, which are the incidence zenith angle, the 

observation zenith angle, and the relative azimuth, respectively. 
The function f(g) in (1) is given by the following ratio: 


Ko(¥o)Ko() 

= 2 
®© Ro(ðo, ð, p) i ? 

where Ro equals to the value of R at zero absorption and 

3 3 7 f 
Ko(¥o0)= i cos Vo + T J dy J Rolo, 0, p) cos? 4 sin Vdd. 
T 
0 0 

(3) 


The function A’o(W) is called the escape function in radiative 
transfer theory. It determines the angular distribution of light es- 
caping the semi-infinite nonabsorbing medium in the framework 
of the Milne problem (with sources located at infinity inside a 
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medium). It can be approximated by the following equation [9] 
at € > 0.2: 


Ko(é) = 5(1 +28) o 
where € = cos Vo. The approximation (4) contributes to biases 
of the ATT as compared to experiments for angles larger than 
approximately 75°. However, the error of (4) is below 2% for 
smaller angles. 

The value of a is given by the following equation: 


lir (5) 


where labs is the absorption path length and l» is the transport 
path length. Studies of the accuracy of (1) as compared to the 
solution of the radiative transfer equation have been conducted 
by Kokhanovsky [6], [7]. It was found that (1) can be applied 
for the values of œ smaller than approximately 1.5. 

Kokhanovsky and Zege [10] derived the following relation- 
ship for the value of œ in the framework of the approximation 
of fractal snow grains: 


a= Byrd (6) 


where y = 47y/A, x is the imaginary part of the refractive 
index of ice, d = 6(V)/(S) is the effective grain size, and 
B x 3.6. Here (V) is the average volume of grains and (S) is 
their average surface area. Clearly, we have for monoidispersed 
spherical grains with the radius a: d = 2a. So d has the meaning 
of a diameter for spheres. 

Combining all equations given above it is possible to arrive 
to the following parameterization of the snow spectral reflection 
function: 


R(90, 9, p) = Ro(o, 9, p) exp(—AV/ yd) (7) 
where 


_ 0.66(1 + cos Jo)(1 + cos ð) 
7 Ro(Wo, 0, p) l 


The constant A is determined by the observation geometry and 
local scattering characteristics of snow. Its dependence on the 
grain size, the wavelength, and the absorption coefficient can be 
neglected. The product c = yd determines the attenuation of the 
optical wave on the grain effective diameter (e.g., it follows for 
the transmitted light intensity: J = To exp(—c), where To is the 
incident light intensity at the surface of the grain). 

Therefore, the snow reflection function is parameterised in 
terms of the diameter d, the reflection function of a semi-infi- 
nite snow layer at zero absorption Ro(Vo, V, y), and the angular 
factor A, which accounts for the influence of geometry on the 
level of light absorption by snow. 

Clearly, the absorption of light by snow should be larger for 
the illumination along normal to a snow layer, and this is ac- 
counted for by the angular factor A. 


A 


(8) 
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Fig. 1. Photo of the measurement site in Hokkaido. 


II. REFLECTION FUNCTION OF SNOW 


A. Experimental Setup 


Measurements of snow reflection functions were made on 
February 9, 2001 (12:00, local time) at Bihoro in eastern 
Hokkaido close to the measurement site described by Aoki 
et al. [2]. The photo of the site is given in Fig. 1. This photo 
was taken not during measurements but at a different day. 

There were no clouds and air temperatures were below 
—10 °C during measurements. The snow depth was 41 cm. The 
snow surface was covered by a thin layer of fresh snow with 
the depth of 2 mm. Most of crystals were dendrites. The snow 
layers to the depth of 7 cm consisted of fine-grained old snow 
and faceted crystals, and the layers from 7 cm to the bottom 
were depth hoar and ice layers. 

Snow grain sizes (radii) were estimated in situ using a hand- 
held lens together with micrographs. This measurement gave 
three kinds of dimensions of grains; one was a half the length 
of the major axis of crystals or dendrites (r4) and the second 
was a half the branch width of dendrites (HBWD) or a half the 
dimension of narrower portion of broken crystals (r2), and the 
third was a snow crystal thickness (r3). We found that r} = 
0.1— 1.0 mm, rə = 0.03—0.15 mm, and r3 = 0.01 — 0.02 mm. 
Aoki et al. [1], [2] concluded that optically equivalent snow 
grain size was rg from the comparison of measured spectral 
albedo with theoretically calculated ones. This means that the 
value of d defined in the previous section must be in the range 
0.06—0.3 mm for the upper layer of snow. 

We find (see below) that the value of r = 0.06 mm (d = 
0.12 mm) must be used to fit measurements using our theory. 
This value falls in the range of measured r2. It confirms the 
finding of Aoki et al. [1], [2] that the optical equivalent snow 
grain size is not rı or r3, but ro (HBWD). 

Our measurements have shown that the snow was contam- 
inated by dust particles and, probably, by other impurities in- 
cluding soot. The contamination was found both at 0—1 cm and 
1-10 cm depths. The newly fallen snow was polluted by dirty 
snow deposited at lower layers. The total mass concentrations 
of impurities was measured and found to be equal to 3.7 ppmw. 
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Fig. 2. Measurement system. 


We did not incorporate the influence of the impurities in the the- 
oretical model because their influence is of importance only in 
visible. Moreover, the refractive index of impurities is not well 
characterized. 

The instruments used for the measurements consist of the 
spectrometer and the pointing system for snow surface. 

The former is a grating spectrometer, “FieldSpec FR,” made 
by ASD Inc. (Boulder, CO), which is the same as that used 
by Aoki et al. [2]. The scanning spectral range of this instru- 
ment was 0.35—2.5 um with the spectral resolution of 3 nm for 
A = 0.35 — 1.0 wm and 10 nm for A = 1.0 — 2.5 um. The 
scanning time is | s with a sampling interval of 1 nm for the full 
spectral range. The detectors are one-dimensional Si photodiode 
CCDs for À = 0.35—1.0 um and two different types of InGaAs 
photodiodes for A = 1.0 — 1.8 um and A = 1.8 — 2.5 um. 

The pointing system is the system to observe the same snow 
surface from any viewing angle (see Fig. 2). The reflected light 
from snow surface is taken into spectrometer by foreoptics with 
the field of view (FOV) equal to 3° and by the optical fiber. 
The foreoptics is connected to the arm which has an angle of 
30° from the rotation axis in the mount, where the rotation axis 
in the mount is set at the same level as snow surface using the 
elevator of tripod. The foreoptics can observe the same snow 
surface from different viewing angles by rotating the arm as 
shown in Fig. 2. The azimuth direction of the bidirectional 
reflectance to be measured is tuneable by changing the azimuth 
direction of the apparatus relative to the sun, where sundial 
is used to detect the solar direction. The downward radiant 
flux, which is necessary to get the reflectance, is measured 
by using white reference standard (WRS) shown in Fig. 2. 
The solar zenith angle was in the range from 53.5° to 54.6° 
during measurements. 

Although the foreoptics points the same snow surface from 
any viewing angle, the precise area of the measured snow sur- 
face varies with the viewing angle. The area of measurements 
is circular for the nadir observations. It changes to an ellipse for 
nonnadir observation geometries. However, the smoother reflec- 
tion functions could be obtained with this system compared with 
the optical system used by Aoki et al. [2]. 
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Fig. 3. Dependence of snow reflection function on the observation zenith 
angle (symbols-measurements for azimuths equal to 0°, 45°, 90°, 180°, 
lines-numerical calculations for the same azimuths at the wavelength 545 
nm. The sun angle was in the range 53.5° to 54.6° during measurements. 
Measurements were taken out at the zenith observation angles equal to 0°, 10°, 
20°, 30°, 40°, 50°, 60°, 70°, 80°, and 85°. 


B. Comparison of Theoretical and Experimental Results 


1) Visible Range: The measured snow reflection function at 
the wavelength 545 nm is given in Fig. 3 for the relative az- 
imuths 0°, 45°, 90°, and 180° (symbols). The azimuth 0° cor- 
responds to the forward scattering in the principal plane. At 
this azimuth, the specular Fresnel reflection from the snow crust 
can appear at solar zenith angles equal to the observation angle. 
There were no crust at the snow surface during measurements 
as it was specified in the previous section. 

Because yd ~ 0 at this wavelength, it follows: R(Jo,0, p) & 
Ro(¥o,0, p). The function Ro(Jo, V, p) has been calculated 
solving the nonlinear integral transport equation (NITE) as 
described by Mishchenko et al. [14]. It was assumed that the 
single-scattering albedo is equal to one and the phase function 
is calculated in the assumption of the fractal snow grains 
(see Macke et al. [12] and Mishchenko et al. [14] for details) 
with the value of d equal to 0.12 mm. Due to the absence of 
absorption and the large size of grains (d >> A), the snow phase 
function does not depend on the value of d (at least outside 
of the small diffraction region, which is of no interest as far 
as light reflection is concerned). Further details on the snow 
particle model used in this paper are given by Macke et al. [12], 
Mishchenko et al. [14], Kokhanovsky [8], and Kokhanovsky 
and Zege [10]. 

The results of numerical calculations for a semi-infinite snow 
layer are given in Fig. 3 by lines. The semi-infinite layer approx- 
imation can be used because of large geometrical thickness of 
snow as described in the previous section. Note that all theoret- 
ical curves were normalized to one value (R = 0.94) obtained 
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from measurements at the observation angle equal to 0° and for 
the solar zenith angle equal to 54°. The theory gives R = 1.0 
in this case, which is by 6% higher than results obtained from 
the experiment for this particular case. This difference could be 
explained by the impurities in snow (e.g., soot and dust). 


Data given in Fig. 3 make it possible to make the following 
conclusions. 


e The radiative transfer equation is capable to predict the 
decrease of the reflection function with the azimuth y as 
found in the experiment. This is mostly due to the decrease 
of the scattering angle 0 and, correspondingly, the increase 
of the phase function with ọ in cases considered. 

° The accuracy is high at azimuths 45° and 90°. 

e The errors of the approximation are small at observation 
angles smaller than 60° for all azimuths studied. 

e The accuracy decreases for the large observation angles 
both for forward scattering (0°) and backward scattering 
(180°). In both cases the radiative transfer equation under- 
estimates the snow reflection. Most probably this is due to 
the specular reflection from the surface of partially oriented 
plate-like crystals positioned in the upper layer of snow. 
The effect is easily observed as a diamond-like shining of a 
newly fallen snow. So the additional component should be 
added to the value of R given by (7). However, such large 
observation angles are not used in snow remote sensing. So 
we do not modify our model to account for possible large 
values of the solar zenith angle. 


Summing up, we conclude that NITE gives an accurate 
approximation for the bidirectional snow reflectance in visible 
(e.g., for the snow satellite remote sensing purposes). Large 
discrepancy of measurements and the theory for the relative 
azimuth 0° and observation angles larger than 60° is not related 
to errors of the AAT as compared to the integro-differential 
radiative transfer equation [14] but rather it is due to conceptual 
problems of the forward model used. 

Indeed, it follows from the experiment that the reflection 
function is close to 1.5 at the observation angle 80° and even 
higher for the observation zenith angle 85° (see Fig. 3). Such 
large values of Fè cannot be described in the framework of the 
standard radiative transfer equation for plane-parallel layers 
[14]. We remind that R gives the ratio of the reflected light 
intensity J to the value of J, if the snow surface is substituted 
by a Lambertian absolutely white surface. Most probably, 
three-dimensional effects and also the specular reflection com- 
ponent must be taken into account to increase the accuracy for 
large zenith observation angles in the principal plane. 

2) Near-Infrared Range: Theoretical curves R(J) ~ Ro() 
given in Fig. 3 have been used [in conjunction with (7) at d = 
0.12 mm] to model the reflection function of snow at several 
wavelengths as is shown in Figs. 4 and 5 for azimuths 45° and 
90°, solar zenith angle 54°, and zenith observation angles 0°, 
10°, 20°, 30°, 40°, 50°, 60°, 70°, 80°, and 85°. 

The value of d is a free parameter of our theory. So it was 
selected to be equal 0.12 mm to have a best fit of our theory 
with respect to measured results. Note that the value of d falls 
in the range of the crystal size variability during measurements 
(see Section III-A). 
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Fig. 4. Dependence of the snow reflection function on the observation zenith 
angle (symbols-measurements at the azimuth equal to 45°, lines-numerical 
calculations) at the wavelengths 545, 1050, 1240, 2210 nm and d = 0.12 mm. 
The values of x have been taken from data given by Warren [17]. The 
illumination/viewing geometry is the same as in Fig. 3. 
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Fig. 5. Same as in Fig. 4 except at the azimuth y = 90°. 


It follows from Figs. 4 and 5 that the exponential factor in 
(7) describes the behavior of experimental curves at azimuths 
45° and 90° for almost all observation angles in the broad spec- 
tral range 545-2210 nm with a high accuracy at d = 0.12 mm. 
The accuracy is poor at observation angles larger than 60° in the 
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Fig. 6. Same as in Fig. 4 except at the azimuth ¢ = 0°. 


principal plane (see Figs. 6 and 7). The principal plane contains 
the Sun, normal to the snow surface and the observation direc- 
tion (the azimuth equals to 0° for forward scattering and the az- 
imuth equals to 180° for backward scattering). This is mostly 
due to the usage of theoretical reflection functions shown in 
Fig. 3. These functions are not accurate enough in the principal 
plane at large observation angles. The accuracy can be consider- 
ably increased (not shown here), if the experimental reflection 
functions shown in Fig. 3 are used in (7) to simulate near-in- 
frared observations. This underlines the fact that the exponen- 
tial factor in (7) captures the physics of the problem in a correct 
way. 

Both the experiment and the theory point out to the fact that 
the dependence of the light reflection on the observation angle 
is stronger for wavelengths, where light absorption is somewhat 
larger. For instance, the slope of F as the function of the obser- 
vation angle in Fig. 4 strongly depends on the level of light ab- 
sorption by snow grains, which is larger for larger wavelengths. 
The accuracy of angular and spectral fits as shown in Fig. 4 is 
quite remarkable taking into account the simplicity of (7) and 
the complexity of light transport in snow. 

Figs. 3—7 confirm that snow cannot be described by a Lam- 
bertian surface model. Then, it follows by definition that R = 
const. This is clearly not the case as indicated in Figs. 3-7. 
The deviations from a Lambertian surface model generally grow 
with the wavelength. 

Minima in experimental curves around the observation angle 
equal 30° in Fig. 5 (and, probably, also maxima around 35° in 
Fig. 4) may be explained by the fact that the snow layer was 
not ideally plane parallel medium as it is assumed in the theo- 
retical model. So these oscillations of experimental curves are 
best explained by the undulation of the snow surface not ac- 
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Fig. 7. Same as in Fig. 4 except at the azimuth y = 180°. 


counted for in the model. The difference in the reflection func- 
tions for azimuths y and ~ + 7 observed in the experiment (not 
shown here) also points in this direction. Actually, the corre- 
spondence of the theory and experiment could be made even 
better if we use as the experimental value not R(y) but R(y) = 
0.5(R(~) + Rl + 7)). Then both experimental errors and the 
inhomogeneity of the surface effects are largely reduced. 

The peaks close to 50° seen in theoretical curves in Fig. 7 are 
not observed in experimental data at y = 180°. These peaks 
originate from the peak of the phase function p(6) of fractal 
particles in the exactly backward direction. This peak is most 
probably due to the inaccuracy of geometrical optics calcula- 
tions for fractal particles at 9 = 180°. So in future for modeling 
snow optical properties the snow phase function obtained in the 
fractal model [12] should be used only for values of 6 < 175°. 
Then one can assume p(@) = p(175°) at 6 € [175°, 180°]. 

The minima in experimental data around 60° (see Fig. 7) are 
due to a shadow caused by the instrument. 


IV. OPTICAL SIZING OF SNOW GRAINS 


We see that (7) can be used to model snow reflection function 
in visible and near-infrared with exclusion of large observation 
angles in the principal plane. This is mostly is due to the fact that 
a is below 1.5 for measurements at wavelengths 545 nm (a = 
0.01), 1050 nm (a & 0.19), 1240 nm (a & 0.43), 2210 nm 
(œ & 1.40) reported here. A high accuracy of the AAT as found 
in this work in a broad spectral range opens up a simple way of 
the determination of the grain size from the measured reflection 
function R by ground-based or orbiting optical instruments. 

In particular, it follows from (7) 


_ both 
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It should be emphasized that all parameters in (9) (except Ro) 
are given by simple analytical formulae specified above. This 
allows for the determination of d from the measured value of 
Rif Ro is obtained from precalculated lookup tables for fractal 
particles. Note that Ro is close to one for fractal particles at the 
nadir observation [10]. If a high accuracy is not of a primary 
importance, then the following parameterization for Ro can be 
used [11]: 


A+ B(E +n) + Cén + p0) 


Ro(Vo, 9,9) = ME +n) (10) 
where 
A=1.247 B=1.186 C=5.157 (11) 
p(0) is the phase function and 
0 = cos™! (-& + /(1— €)(1 — n?) cos 9) (12) 


is the scattering angle, and € = cos Vo, n = cos V. The calcu- 
lation of p(@) for irregularly shaped ice crystals requires long 
computation time. Therefore, the following parameterization 
has been developed [11]: 

p(0) = 11.1 exp(—0.0876) + 1.1 exp(—0.0146) (13) 
where 6 > 10° is expressed in degrees. The accuracy of (10) 
is better than 3% for the nadir observation and incidence angles 
smaller than 78°. The error is below 10% at € > 0.2, n > 0.9. 

If the value of d is obtained from measurements in the near 
infrared (e.g., at 1240 nm), then measurements in visible (e.g., at 
545 nm) can be used to deduce the snow absorption coefficient 
vis. Namely, it follows from (9) 


1 R 
= — hR’ ( — }. 
Ad Ro 
Clearly, Yvis should correlate with the soot concentration c (and, 
possibly, with abundances of other absorbers) as discussed by 


Zege et al. [15]. This allows for the development of techniques 
to monitor the contamination of snow fields by soot from space. 


Yois (14) 


V. CONCLUSION 


The asymptotic radiative transfer theory as applied to the 
snow reflection function calculation [10] has been used for the 
interpretation of in situ measurements performed in Hokkaido, 
Japan on February 9, 2001. A close correspondence of theory 
and experimental results is found outside of the principal plane 
using the effective grain size d = 0.12 mm. This means that 
we also retrieved the optically equivalent size of snow grains 
by comparing our measurements with the theory. The value of 
d falls in the range of measured lengths of the branch width of 
dendrites or dimensions of narrower portions of broken crys- 
tals. This confirms similar findings of Aoki et al. [1], [2]. More 
measurements are needed to establish the correlation between 
the optically equivalent radius as retrieved from the reflectance 
measurements with actual sizes of grains as measured in situ. 
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In the principal plane, which holds the observation and illu- 
mination directions, the accuracy of the approximation is low 
for large observation angles. This can be attributed to effects of 
the specular light reflections by tiny ice crystals deposited on the 
snow surface. These reflections are not accounted in the model 
considered here. 

We underline that only a small subset of measurements per- 
formed was analyzed in this work. However, we found that the 
subset used captures the major angular characteristics of the 
snow reflectance. 

Results obtained show that our theory [see (7)] can be used 
to monitor the grain size and also the level of snow pollution 
from spaceborne optical instrumentation. For this, however, the 
atmospheric correction of snow imagery should be performed, 
which is not an easy task taking into account problems with the 
aerosol optical thickness retrieval over bright surfaces. How- 
ever, this is a crucial point only in visible. The snow is almost 
black in near-infrared (see Figs. 4-7). This makes retrieval of the 
aerosol properties over snow in near-infrared much more simple 
than this is the case in the visible. 

Other urgent point to be addressed is the relation of the 
retrieved optically equivalent ice grain size to that measured 
in situ. For this, however, dedicated campaigns with much 
more advanced techniques as compared to those used in this 
study for in situ ice grains particle sizing must be utilized (e.g., 
in conjunction with spectral reflectance measurements). In 
particular, the automatic image analysis techniques look very 
promising in this respect. 
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